********************************************************************************
********************************************************************************
* Title: 		cortisol
* Project:		Stress, Ethnicity, and Prosocial Behavior
* Author:		Moritz Poll (moritz.poll@brown.edu)
* PIs:			Johannes Haushofer, Sara Lowes, Abednego Musau, David Ndetei, 
*				Nathan Nunn, Moritz Poll, Nancy Qian
* Purpose:		Cortisol Analysis Visualization
********************************************************************************
********************************************************************************

set more off
noisily dis as input "Start: Cortisol analysis"

********************************************************************************
***************************** Average game timing ******************************
********************************************************************************

use  "$cleandata_dir/HIO_for_analysis_session.dta", clear
include "$do_analysis_dir/definitions.do"

local sformat `" "%9.0f" "'
* Start and end
sum ts_task1_start_pill if ts_taskorder == 1
local  dg1s = string(r(mean), `sformat')
sum ts_task1_end_pill   if ts_taskorder == 1
local  dg1e = string(r(mean), `sformat')
sum ts_task2_start_pill if ts_taskorder == 1
local cyd1s = string(r(mean), `sformat')
sum ts_task2_end_pill   if ts_taskorder == 1
local cyd1e = string(r(mean), `sformat')
sum ts_task3_start_pill if ts_taskorder == 1
local tg11s = string(r(mean), `sformat')
sum ts_task3_end_pill   if ts_taskorder == 1
local tg11e = string(r(mean), `sformat')
sum ts_task4_start_pill if ts_taskorder == 1
local tg21s = string(r(mean), `sformat')
sum ts_task4_end_pill   if ts_taskorder == 1
local tg21e = string(r(mean), `sformat')
sum ts_task1_start_pill if ts_taskorder == 2
local tg12s = string(r(mean), `sformat')
sum ts_task1_end_pill   if ts_taskorder == 2
local tg12e = string(r(mean), `sformat')
sum ts_task2_start_pill if ts_taskorder == 2
local tg22s = string(r(mean), `sformat')
sum ts_task2_end_pill   if ts_taskorder == 2
local tg22e = string(r(mean), `sformat')
sum ts_task3_start_pill if ts_taskorder == 2
local  dg2s = string(r(mean), `sformat')
sum ts_task3_end_pill   if ts_taskorder == 2
local  dg2e = string(r(mean), `sformat')
sum ts_task4_start_pill if ts_taskorder == 2
local cyd2s = string(r(mean), `sformat')
sum ts_task4_end_pill   if ts_taskorder == 2
local cyd2e = string(r(mean), `sformat')
sum ts_exit_start_pill
local social_proximity_s = string(r(mean), `sformat')
sum ts_exit_end_pill
local social_proximity_e = string(r(mean), `sformat')
sum ts_lab_start_pill
local instructions_s = string(r(mean), `sformat')
sum ts_salivette_pill1
local instructions_e = string(r(mean), `sformat')
sum ts_salivette_pill1
local pause_s = string(r(mean), `sformat')
sum ts_salivette_pill2
local pause_e = string(r(mean), `sformat')

* Duration
foreach game in pause instructions $main_table_vars {
	forvalues order = 1 / 2 {
		capture local `game'`order'd = ``game'`order'e' - ``game'`order's'
		capture local `game'_d = ``game'_e' - ``game'_s' // for social proximity I disregard order
	}
}

********************************************************************************
********************** Graph over salivette/time with SEs **********************
********************************************************************************

* Plot cortisol levels against salivettes
use "$cleandata_dir/HIO_for_analysis_wide.dta", clear

forvalues s = 1 / 6 {
	ttest cortisol`s', by(treatment)
	local s`s' = "n.s."
	if r(p) < 0.1  local s`s' = "*"
	if r(p) < 0.05 local s`s' = "**"
	if r(p) < 0.01 local s`s' = "***"
	
	gen cortisol_r`s' = cortisol`s'/cortisol1
	label var cortisol_r`s' "Cortisol level at salivette `s' relative to individual baseline"
	
	ttest cortisol_r`s', by(treatment)
	local sr`s' = "n.s."
	if r(p) < 0.1  local sr`s' = "*"
	if r(p) < 0.05 local sr`s' = "**"
	if r(p) < 0.01 local sr`s' = "***"
	
	sum cortisol1
	gen cortisol_z`s' = (cortisol`s' - r(mean))/r(sd)
	label var cortisol_z`s' "Cortisol level (z-scored)"
	
	ttest cortisol_z`s', by(treatment)
	local sz`s' = "n.s."
	if r(p) < 0.1  local sz`s' = "*"
	if r(p) < 0.05 local sz`s' = "**"
	if r(p) < 0.01 local sz`s' = "***"
}


collapse cortisol* salivettetime* ts_taskorder (semean) cse1=cortisol1 cse2=cortisol2 cse3=cortisol3 cse4=cortisol4 cse5=cortisol5 cse6=cortisol6 cser1=cortisol_r1 cser2=cortisol_r2 cser3=cortisol_r3 cser4=cortisol_r4 cser5=cortisol_r5 cser6=cortisol_r6 csez1=cortisol_z1 csez2=cortisol_z2 csez3=cortisol_z3 csez4=cortisol_z4 csez5=cortisol_z5 csez6=cortisol_z6, by(treatment)
reshape long cortisol cortisol_r cortisol_z cse cser csez salivettetime salivettetime_pill, i(treatment) j(salivette)

gen csel = cortisol - 1.96*cse
gen cseh = cortisol + 1.96*cse
gen cselr = cortisol_r - 1.96*cser
gen csehr = cortisol_r + 1.96*cser
gen cselz = cortisol_z - 1.96*csez
gen csehz = cortisol_z + 1.96*csez

sum cseh
local yhigh = r(max)
sum csehr
local yhighr = r(max)
sum csehz
local yhighz = r(max)
		
local colorP "edkblue"
local colorH "maroon"
twoway  (line cortisol salivette if treatment == 0, color(`colorP') lpattern(dash)) (rcap csel cseh salivette if treatment == 0, color(`colorP')) ///
		(line cortisol salivette if treatment == 1, color(`colorH')) 				(rcap csel cseh salivette if treatment == 1, color(`colorH')), ///
		text(210 2.75 "                                                                                                                   ", box bcolor(gs14)) ///
		text(200 2.75 "                                                                                                                   ", box bcolor(gs14)) ///
		text(220 2.75 "                                                                                                                   ", box bcolor(gs14)) ///
		text(230 2.75 "                                                                                                                   ", box bcolor(gs14)) ///
		text(240 2.75 "                                                                                                                   ", box bcolor(gs14)) ///
		text(232.5 -0.75 "Task order 1:", color(gs0) box bcolor(none)) text(207.5 -0.75 "Task order 2:", color(gs5) box bcolor(none)) text(225 0.6 "Instructions", color(gs3) box bcolor(none)) text(225 1.6 "Break", color(gs3) box bcolor(none)) ///
		text(215 0.6 "`instructions_d' min", color(gs3) box bcolor(none)) text(215 1.6 "`pause_d' min", color(gs3) box bcolor(none)) ///
		text(237.5 2.5 "DG" , color(gs0) box bcolor(none)) text(237.5 3.5 "CYD", color(gs0) box bcolor(none)) text(237.5 4.5 "TG1", color(gs0) box bcolor(none)) text(237.5 5.5 "TG2", color(gs0) box bcolor(none)) ///
		text(212.5 2.5 "TG1", color(gs5) box bcolor(none)) text(212.5 3.5 "TG2", color(gs5) box bcolor(none)) text(212.5 4.5 "DG" , color(gs5) box bcolor(none)) text(212.5 5.5 "CYD", color(gs5) box bcolor(none)) ///
		text(227.5 2.5 "`dg1d' min" , color(gs0) box bcolor(none)) text(227.5 3.5 "`cyd1d' min", color(gs0) box bcolor(none)) text(227.5 4.5 "`tg11d' min", color(gs0) box bcolor(none)) text(227.5 5.5 "`tg21d' min", color(gs0) box bcolor(none)) ///
		text(202.5 2.5 "`tg12d' min", color(gs5) box bcolor(none)) text(202.5 3.5 "`tg22d' min", color(gs5) box bcolor(none)) text(202.5 4.5 "`dg2d' min" , color(gs5) box bcolor(none)) text(202.5 5.5 "`cyd2d' min", color(gs5) box bcolor(none)) ///
		text(230 6.5 "Social", color(gs3) box bcolor(none)) text(220 6.5 "Proximity", color(gs3) box bcolor(none)) text(210 6.5 "`social_proximity_d' min", color(gs3) box bcolor(none)) ///
 		text(226 4 "_______________________", color(white)) text(226 -0.75 "________", color(white)) ///
		xline(1.25, lcolor(gs10) lpattern(shortdash) lwidth(medium)) text(`=`yhigh'+10' 1 "`s1'", color(gs8)) text(`=`yhigh'+10' 2 "`s2'", color(gs8)) text(`=`yhigh'+10' 3 "`s3'", color(gs8)) text(`=`yhigh'+10' 4 "`s4'", color(gs8)) text(`=`yhigh'+10' 5 "`s5'", color(gs8)) text(`=`yhigh'+10' 6 "`s6'", color(gs8)) ///
		xscale(r(0 6.75)) xlabel(1(1)6 1.25 `""|" "Pill""') xtitle(Salivette Number) ylabel(0(50)150) yscale(r(0 180)) ytitle(Salivary cortisol in nmol/l) legend(order(1 "Placebo" 3 "Hydrocortisone") symxsize(12) region(lcolor(white))) graphregion(color(white) margin(12 3 3 17))
graph export "$figures_dir/cortisol_salivette.pdf", replace

noisily dis as input "End: Cortisol analysis"
